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ABSTRACT 

The equations of secular evolution for dust grains in mean motion resonances with 
a planet are solved for stationary points. This is done including both Poynting- 
Robertson effect and stellar wind. The solutions are stationary in semimajor axis, 
eccentricity, and resonant angle, but allow the pericentre to advance. The semimajor 
axis of stationary solutions can be slightly shifted from the exact resonant value. The 
periodicity of the stationary solutions in a reference frame orbiting with the planet 
is analytically proved. The existence of periodic solutions in mean motion resonances 
means that analytical theory enables for dust particles also infinitely long capture 
times. The stationary solutions are periodic motions to which the eccentricity asymp¬ 
totically approaches and around which the libration occurs. 

Using numerical integration of equation of motion are successfully found initial 
conditions corresponding to the stationary solutions. Numerically and analytically 
determined shifts of the semimajor axis form the exact resonance for the stationary 
solutions are in excellent agreement. 

The stationary solutions can be plotted by locations of pericenters in the reference 
frame orbiting with the planet. The pericenters are distributed in the space according 
to properties of dust particles. 

Key words: interplanetary medium - zodiacal dust - stars: general - celestial me¬ 
chanics 


1 INTRODUCTION 

A lenticular cloud of dust centered on the Sun with its main 
axis lying in the ecliptic plane creates well known zodiacal 
light. Discoveries and improving observations of debris disks 
around main-sequence stars (e.g. HoUande£aL|1998[^ Debes 


Weinberger & Kuchner 2009; Schneider et al. 2014; Lagrange 
et al. 12016 1 motivate to understand better the dynamics of 

the dust particles in such environments. Sharp inner edges, 
clumps and asymmetries observed in the debris disks are 
commonly interpreted as signatures of a planet created in 
earlier stages of the disk evolution and located near or within 
the disks (e.g. Wyatt|2003 Kuchner fe Holman|2003 Stark 


fe Kuchner||2009~ Rodigas, Malhotra fe Hinz||2014 l. Such a 

planet can gravitationally capture the dust particles in the 
so called mean motion resonances (MMRs) and change the 
morphology of the disk. 

* pavol.pastor@hvezdarenlevice.sk 


Mean motion resonances occur when a ratio of orbital 
periods of the dust particle and the planet is near a ra¬ 
tio of two small natural numbers. The MMRs in a circu¬ 
lar restricted three-body problem (CRTBP) with a stellar 
electromagnetic radiation and a stellar wind include numer¬ 
ous aspects which are directly applicable to astrophysical 
systems, ranging from the resonant dust ring close to the 


Earth’s orbit (e.g. Jackson & Zook 1989 |Dermott et al. 
1994 Reach et al.|1995 Reach|2010 Pastor|2014b I through 


dust in the Edgeworth-Kuiper belt zone (e.g. Liou & Zook 
1999| |Moro-Martm fe Malhotra||2002[ |Holmes et al.||2003[ 
Kuchner & Stark 20101 to structures observed in debris 


disks (e.g. Kuchner fe Holman|2003 Stark fe Kuchner|2008 


Shannon, Mustill & Wyatt 20151. Analytical results can be 
applied into semi-analytic models of resonant dust just to 
improve computation speed of resonant disk patterns (e.g. 
Wyatt||2003| |Kuchner fe Holman] |2003| Stark fe Kuchner| 


2009| Rodigas et al. ||20 1 -T[~ Shannon et al. 20151. 

In this work we show that the CRTBP with radiation 
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2 P. Pastor 


comprises also periodic solutions for the motion of the dust 
particle in an exterior MMR (the orbital period of the dust 
particle is larger than the orbital period of the planet in the 
exterior MMRs). The periodic solutions exist at “univer¬ 
sal eccentricities” found by Beauge & Ferraz-Mello 1994 1 . 


The eccentricity of dust orbit during a capture in the exte¬ 


rior MMR approaches the universal eccentricity (e.g. Liou 


fe Zook|1997 Klacka et al.|2008 Stark fc Kuchner|2009 1. If 

such periodic solutions exist, then maximal capture times for 
the dust particles in the exterior MMRs should be infinite. 


Using adiabatic invariant theory Gomes (19951 found 


a criterion determining whether the libration amplitude in¬ 
creases or decreases for dissipative forces which lead to an 
evolution of the eccentricity characterised with asymptotic 
approach of the eccentricity to the universal eccentricity. 
Results of averaged theory in |Gomes (1995]) showed that 
the libration amplitude must increase for the dust particles 
in all exterior resonances in the planar CRTBP with the 
Poynting-Robertson (PR) effect. |Gomes| (1995] results al¬ 
low also the existence of the periodic solutions. 

Periodic solutions in 1/1 resonance were already found 
to exist ( |Liou, Zook fe Jackson||l995| |Pastor|[2014b[ ). The 
existence of secularly stationary solutions in the vicinity 
of the triangular Lagrangian equilibrium points in spatial- 
circular, planar-elliptic and spatial-elliptic restricted three- 
body problem with the PR effect was recently confirmed in 
|Lhotka fe CellettI (20151. 


The content of this paper is the following. In Section 
[2] are presented averaged resonant equations for the MMRs 
in the planar CRTBP with arbitrary non-gravitational ef¬ 
fects for which secular theory can be applied (general case). 
In Section [ 3 ] are determined conditions for the stationary 
points in the general case. In Section [I] is shown that for 
astrophysical systems with rotational symmetry an evolu¬ 
tion during the stationary solution in the MMR is always 
periodic. In Section[5] are derived the conditions for the sta¬ 
tionary points when the non-gravitational effects are the PR 
effect and the radial stellar wind. In Section [6] we determine 
the stationary solutions for the considered non-gravitational 
effects, we discuss properties of the solutions and present 
pericenters of the orbits corresponding to the stationary so¬ 
lutions in a reference frame orbiting with the Earth. Con¬ 
clusions are drawn in Section [3 


2 AVERAGED RESONANT EQUATIONS FOR 
DUST PARTICLE INFLUENCED BY 
NON-GRAVITATIONAL EFFECTS 

Orbital evolution in a dynamical system is often averaged 
over a period of the dynamical system (secular framework). 
The secular framework also enables define a stationarity. In 
the secular framework the stationarity is obtained if one or 
more orbital elements remain constant after the averaging 
process. In the MMRs the semimajor axis of the particle 
orbit is also secularly stationary after the averaging over an 
libration period. In what follows the word secularly before 


the stationary will be omitted. A dust particle can be cap¬ 
tured in an MMR with a planet in the planar CRTBP and 
simultaneously be affected by non-gravitational effects. In 
the general case secular variations of a particle orbit, caused 
by the non-gravitational effects, can depend on spatial orien¬ 
tation of the orbit. In this case resonant equations averaged 
over a synodic period have form (Pastor|[2014a): 
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Here, k = ecos a and h = e sin a are the non-canonical res¬ 
onant variables, e is the eccentricity of the particle orbit, 
and a is the resonant angular variable cr = ip/q — Co. The 
expression for the resonant angular variable contains: ip = 
(p + q) Ap — pX, two integers p and q, the mean longitude A, 
and the longitude of pericenter Co. The subscript P is used 
for quantities belonging to the planet, a = \/l~ e 2 , s = 
p/q, and <ji = ip/q — Cop. R is the disturbing function of 
the planar CRTBP. L = y/pa, where, p = GqM*, Gq is the 
gravitational constant, M* is the mass of the star, and a is 
the semimajor axis of the particle orbit, n = yf p/a 3 is the 
mean motion of the particle. d*a denotes the partial deriva¬ 
tive with respect to the semimajor axis calculated with an 
assumption that the mean motion of the particle n is not 
a function of the semimajor axis (see e.g. jPanb y 1 988] . We 
use angle <7b defined so that the mean anomaly M can be 
computed from M = nt + <Jb (Bate, Mueller & Wh ite|l97 1). 
In this case the following relations hold between the aver¬ 


aged partial derivatives of the disturbing function (Pastor 


2014a) 
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From the definitions after Eqs. |l] we have a = ip/q — Co 
and ay = ip/q — Cop. Therefore, the last equation in Eqs. 0 
reads 


da i 
dt 


da duo 
dt dt 


( 3 ) 


where we have used the fact that Cop = 0 in the CRTBP. 
When we again use the definitions after Eqs. 0 for L = 
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Locations of stationary/periodic solutions 3 


y/jla, k = e cos a and h = e sin a we can obtain, using some 
analytic manipulations, the following useful relations 
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Eq. @ and Eqs. 0 can be used in order to rewrite 
as follows 
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The system given by Eqs. 0 enables to study the secular 
orbital evolution of the dust particle captured in the specific 
MMR given by the resonant numbers p and q under the ac¬ 
tion of the non-gravitational effects. The expression of Eqs. 
0 using the orbital elements and the resonant angular vari¬ 
able, given by Eqs. 0, will be useful for a determination of 
properties of the stationary solutions as we will see later. 
For the averaged non-conservative terms (Qi) and (Pi), i = 
1, 2 hold ( Pastor|[2014a l 


«w~<£>_-<£ + *£>, 

«»-(*}.■ 


{P 2 } = - 
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dt / BF 


na 2 e / de \ 
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( 6 ) 


where (da/dt) ef, (de/dt) ef, (dui/dt) ef, and (dah/dt + t 
dn/dt)EF are caused by the non-gravitational effects only. 
The averaged values on the right-hand sides in Eqs. § can 
be obtained using the Gaussian perturbation equations of 
celestial mechanics (e.g. Danby 19881. Their specific values 
for the PR effect and the radial stellar wind will be shown 
later. The equations in this section are derived for the gen¬ 
eral case. 


3 CONDITIONS FOR STATIONARY POINTS 


If Eqs. 0 hold in a case without the non-gravitational ef¬ 
fects (conservative problem), then they affect also the evo¬ 
lution of the semimajor axis. This is due to the fact that the 
condition de/dt = 0 implies da/dt = 0. This can be most 
easily seen if we calculate time derivative of the eccentricity 
averaged over the synodic period using the first two equa¬ 
tions in Eqs. 0 and the last two equations in Eqs. 0 (see 
also Eq. 17 in Pastor|[2014a) 
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In the conservative problem Eq. 0 gives for de/dt = 0 
= 0 or dL/dt = 0 (see also |Beauge fe Ferraz-| 
I. Eq. 0 shows that in the case with the non- 
gravitational effects de/dt = 0 is not sufficient for dL/dt = 
0, in general. We are interested in the stationary solutions of 
the system of equations given by Eqs. 0 restricted by the 
same conditions as in the conservative problem. If we add 
the condition dL/dt = 0 to the conditions in Eqs. 0, then 
the semimajor axis will be also stationary similarly as in 
the conservative problem at the stationary points in the kh 
plane. The condition dL/dt = 0 is necessary for the period¬ 
icity of the stationary solutions as we will see later. Hence, 
the stationary points will be determined by the substitution 
of the following three equations into Eqs. 0 


that da/dt 


Mello 1994 


dk dh dL 
dt dt dt 


(9) 


Eqs. (|9j) were used also in 

Beauge & Ferraz-Mello (1993|) 

Sidlichovsky & Nesvorny 

(19941. 


The substitution of Eqs. 0 into Eqs. 0 leads to the 
conditions 
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as can be seen from Eqs. 0. The semimajor axis and the 
position of the libration center in the kh plane are stationary 
in the sought for stationary points. 


4 PERIODICITY OF EVOLUTIONS 

CORRESPONDING TO STATIONARY 
SOLUTIONS IN RESONANCES 


The stationary points in the kh plane are determined by 


dk dh 
dt dt 


(7) 


In the exact resonance the averaged mean motion of the par¬ 
ticle and the mean motion of the planet have a ratio that 
is exactly equal to the ratio of two natural numbers. The 
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4 P. Pastor 


periodicity of evolutions corresponding to the stationary so¬ 
lutions in the exact resonances with fixed longitude of peri- 


center was already showed in Pastor (2014a). In this work we 
take into account the possibility that the longitude of peri- 
center can be variable. Influence of the variable longitude 
of pericenter on the periodicity of the evolutions during the 
stationary solutions will be determined. Beside the variable 
longitude of pericenter in this work we consider also non¬ 
exact resonances. In this section, the symbols used denote 
the un-averaged quantities. We assume that the physical pa¬ 
rameters of the planar CRTBP with the non-gravitational 
effects do not change during the synodic period. In our case, 
the un-averaged resonant angular variable is 


P + 9, > - 

a = -Ap — sA — lo 


with the mean longitude 


A — Af -f- to — nt -)— o" p, —(— cx_? . 


( 12 ) 


(13) 


For the time derivative of the mean longitude in Eqs. (13I 
we obtain 


dX dab dn dio 

dt dt dt dt 
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Using the previous equation we obtain for the time deriva¬ 
tive of a in Eq. (12I 
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Let the motion be stationary according to the conditions in 


Eqs. (111. Eq. (151 can be averaged over the synodic period. 


If we use the last condition in Eqs. (Ill for the stationarity 


of the resonant angular variable in the result of averaging, 
we obtain 
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In |Pastor| ( |2014a| ) only the exact resonances were con¬ 
sidered. The condition np [p + q)/p = (1/Ts) /„ Ts n dt valid 
for the exact resonances may not be valid for the non-exact 


resonances. Therefore, the first two terms in Eq. (161 may 


not cancel each other as in the case for the exact resonances. 

The periodicity of evolution corresponding to the sta¬ 
tionary solution in P astor| ( |201 4a| was demostrated for the 
non-gravitational effects which include an interstellar wind. 
For the star moving through a cloud of an interstellar mat¬ 
ter the interstellar wind is monodirectional (for details see 
Baines, Williams & Asebiomo 1965). The monodirectional 


character of the interstellar wind enables only the stationary 
solutions for which the longitude of pericenter is stationary. 


If we substitute Eq. (161 in the averaged Eq. (141, then we 
obtain 


1 f Ts dX , p + q q 1 f Ts dio 
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Ts J 0 dt p p T s J 0 dt 


(17) 


In the equation above we can see that the difference between 
mean longitudes after the synodic period Ts depends on the 
synodic period and on the averaged time derivative of the 
longitude of pericenter. In what follows we determine a rela¬ 
tion between these two variables. Eqs. 0 were obtained by 
the averaging over the synodic period using a cyclic angular 
variable ap as the integration v ariable, ao is defined (see e.g. 


Sidlichovsky fe Nesvorn/p9941 as follows 
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The difference between 02 at time zero and 02 after one 
synodic period Ts is equal to 27r. Implicit differentiation of 


<T 2 defined in Eq. (181 yields 
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The last equation can be integrated over the synodic period 


2n= - 

q [Jo 


J 

Uo 

f 

Jo 


n dt + 


r T s 

Jo 


(dab dn\ , 
i —+ + t — )dt 
\ dt dt J 


t Ts did 

— dt - n P T s 
at 


( 20 ) 


Substitution Eq. (161 into Eq. (20) yields the sought for 


relation between the synodic period and the averaged time 
derivative of the longitude of pericenter 
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The difference between mean anomalies after one synodic 
period can be obtained using 


kCXL 
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— dt + — 
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and Eq. (21 1 in Eq. 1171. We obtain 


1 f ls dM , , . 2tt 

%L ur dt = {r +- ) t. 

The last equation finally implies 

-M(Ts) — M(0) = (p + q) 2n 


Now, we rewrite Eq. (21 1 as follows 


<I>(0) — Ap(0) = w(Ts) — Ap(Ts) + 2np 


( 22 ) 


(23) 


(24) 


(25) 


where we have used the identity Ap(Ts) — Ap(0) = npTs- 
Eq. ( |25[ ) shows that a difference in the angle between the 
planet and the pericenter of the particle orbit after the syn¬ 


odic period is an integer multiplied by 2-7T. Eq. (251 takes into 


The third term in Eq. (161 is caused by the variability (in- account the arbitrary variability of the longitude of pericen- 


stationarity) of the longitude of pericenter and in Pastor 
( |2014aD was equal to zero. 

Eq. (|14|) can also be averaged over the synodic period. 


ter. Eq. (241 gives that after the synodic period the particle 


is in the same position on the orbit as in the time zero. 
The first two equations in Eqs. © ensure that after the 
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Figure 1 . Evolution of a dust particle with R ^ = 5 ^xm, g = 2 g/cm 3 , and Qp r = 1 captured in an exterior mean motion 6/5 orbital 
resonance in the planar CRTBP with the PR effect and the radial solar wind which comprises the Earth and the Sun as the two main 
bodies. Depicted evolving parameters are described in the text. The dashed line in the second plot shows the asymptotic value of the 
eccentricity. The dashed line in the kh plane indicates the values of e and a that lead to collisions of the planet with the particle. All 
evolving parameters are averaged over the synodic period. 


Table 1. The universal eccentricity (the first equation in Eqs. |3l[ ) for several exterior resonances. The eccentricities implicate from Eqs. 
when the non-gravitational effects are the PR effect and the radial stellar wind. 
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A = 0 AU 


A = 3 x 10 -6 AU 
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Figure 2. Solutions of the conditions for the stationary points in the planar CRTBP with radiation (Eqs. |3l| ). Solutions of the second 
equation, determined at the universal eccentricity, are depicted with red curves. Similarly solutions of the third equation are depicted 
with green curves. Dust particles are captured in an exterior 6/5 MMR with the Earth in a circular orbit around the Sun. In the left 
plot is A = 0 used and the right plot was obtained using A = 3 x 10 —6 AU. Black circles indicate intersections of red and green lines 
representing solutions of the system given by Eqs. ©■ 


synodic period the semimajor axis and the eccentricity are 
the same as in the time zero. On the basis of Eqs. GD- 
Eq. (241 and Eq. (251 we can draw the following conclusions 
for the evolution during the stationary solution in the reso¬ 
nance. If the secular time derivatives of the orbital elements 
caused by the non-gravitational effects (non-gravitational 
secular variations) do not depend on the longitude of peri- 
center, then the evolution will be repeated after the synodic 
period (the motion of the dust particle is periodic). The 
non-gravitational secular variations do not depend on the 
longitude of pericenter in astrophysical systems with rota¬ 
tional symmetry (Appendix |A| . In the cases when the non- 
gravitational secular variations depend on the longitude of 
pericenter the variation of the longitude of pericenter after 
one synodic period can cancel the periodicity. 


5 CONDITIONS FOR STATIONARY POINTS 
WITH PR EFFECT AND RADIAL STELLAR 
WIND 


Let us consider a star which produces an electromagnetic ra¬ 
diation and a radial stellar wind with one planet in a circular 
orbit. In astrophysical modelling of the orbital evolution of 
dust grains it is often assumed that a spherical particle can 
be used as an approximation to real, the arbitrarily shaped, 
particle. The force associated with the action of the electro¬ 
magnetic radiation on the moving the spherical dust parti¬ 
cle is the PR effect (|Poynting|1904| |Robertson|1937| |Klacka| 


2004 Klacka et al.|[2014 |. The stellar wind also affects the 

dynamics of the dust particles. A relativistically covariant 
equation of motion of the dust particle caused by a parallel 


beam of corpuscules impinging on the dust particle was de¬ 
rived in |Klack a et al. (2012]). On the basis of the mentioned 
papers we obtain the following equation of motion for the 
spherical dust particle in a reference frame associated with 
the star in accuracy to first order in v/c (v is the speed of 
the dust particle with respect to the star and c is the speed 
of light in vacuum), first order in u/c (u is the speed of the 
stellar wind with respect to the star) and first order in v/u 


dv 

dt 


-^(l-/?)e R - 


GqMp , _ _ , GoMp _ 

-3 ( r " rp ) “ rp 




1 + 


<4 


r — rp 
v ■ e R 


v 

eR H — 
c 


(26) 


where r is the radial distance between the star and the dust 
particle, e*R is the unit vector directed from the star to the 
particle, v the velocity of the particle with respect to the 
star, Mp is the mass of the planet, rp is the position vec¬ 
tor of the planet with respect to the star, and rp = |rp|. 
The parameter (3 is defined as the ratio between the elec¬ 
tromagnetic radiation pressure force and the gravitational 
force between the star and the particle at rest with respect 
to the star. For /3 we have from its definition 


3L*Qpr 

16ncfiRdQ ’ 


(27) 


where L* is the stellar luminosity, Q pl is the dimensionless 
efficiency factor for the radiation pressure averaged over the 
stellar spectrum and calculated for the radial direction ( Q pl 
= 1 for a perfectly absorbing sphere), c is the speed of light 
in vacuum, and Rd is the radius of the dust particle with 
mass density g. r/ is to the given accuracy the ratio of the 
stellar wind energy to the stellar electromagnetic radiation 
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Figure 3. Evo lution of orbital parameters averaged over the synodic period obtained from numerical integration of the equation of 
motion (Eq. |26[ > close to the stationary point 12b depicted in Fig. |3J A dust particle with S = 0.02027611 and Q pI = 1 is captured in an 
exterior mean motion 6/5 orbital resonance with the Earth in a circular orbit around the Sun under the action of the PR effect and the 
radial solar wind. The stationary solution occurs in a non-exact MMR. Un-averaged evolution during the first 30 years can be seen in 
Fig. [4] 


energy, both radiated per unit time 
. 2 N 

47rr u 2 

V / j ^sw i^Tl sw ic 
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(28) 


where m s 


and n s . 


i = 1 to N, are the masses and 


concentrations of the stellar wind particles at a distance r 
from the star (u = 450 km/s and r/ = 0.38 for the Sun, 
jKlacka et al.|2012fr . 

The orbital evolution will be given by the evolution of 
osculating orbital elements calculated for the case when a 
cent ral a cceleration is defined by the first Keplerian term in 
Eq. (26 1, namely —/z(l — /3) en/r 2 . The expressions in the 
previous sections which contain L and n must be modified 
for this central acceleration using a reduced mass of the star 
M*( 1 — P). The modifying equations are L = \J /r(l — /3 )a 
and n = \Jn(l — fi)/a 3 . Using the last term in Eq. (26 ) 


as a perturbation of the orbital motion we obtain from the 
Gaussian perturbation equations of celestial mechanics the 
following averaged values in Eqs. © 
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The substitution of Eqs. (291 into Eqs. © gives 

(Qi) = o , 

(Q2) = 0 , 


(Pi) = - 


Pfj,n 
2 ca 3 


1 + J”) ( 3e2 + 2) 

Vpr/ 


(P>) = -(Pi)- P -f 


1 + S- 
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(30) 


Now, we return to the conditions for the stationary 
points during a capture in the MMRs occurring in the planar 
CRTBP with non-gravitational effects given by Eqs. ITo|. 


The substitution of Eqs. (301 in Eqs. (101 yields for the PR 


effect and the radial solar wind 
3e 2 + 2 p + q 
2 a 3 p 
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— 0.24724 
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Figure 4. Evolution of osculating orbital parameters during initial 30 years of the evolution depicted in Fig. [ 3 ] (black line). For the sake 
of comparison also the evolution of orbital parameters averaged over the synodic period is shown (blue line). The averaged evolution is 
created only by four points placed in the middle of a given averaged synodic period. Scales on vertical axes are different from the scales 
in Fig.[3]in order to show full evolution of osculating a, e and a during the synodic period. 


a r , dR p + q 2 sa OR 

— 1 + s (1 - a) Hr- = -—-n P -sn+— — 
Le oe q L o*a 


dR 

da 


q ppn 
p + q c 


1 + PL 

Q' P r 


(31) 


For an interior resonance the first equation in Eqs. (311 does 


not have any solution, p + q < p holds in an exterior reso¬ 


nance and the solution of the first equation in Eqs. (311 is 
the “universal eccentricity” (|Beauge fe Ferraz-Mcllo 


1994 


Liou fe Zook|[l997|. The universal eccentricities for several 


exterior resonances are in Table |T] 

An example of the evolution approaching to the uni¬ 
versal eccentricity is depicted in Fig. [I] The evolution is 
obtained from numerical solution of the equation of motion 
(Eq. 26 1 for a dust particle with Rd = 5 /im, q = 2 g/cm 3 , 
and <2 pr = 1 captured in an exterior mean motion 6/5 orbital 
resonance with the Earth in a circular orbit around the Sun 
under the action the PR effect and the radial solar wind. 
The depicted parameters are from left to right and from top 
to bottom, as follows: a (semimajor axis), e (eccentricity), 
uj (longitude of perihelion), a (resonant angular variable), 
kh point, k = ecos a, h = esin a, L = yj /r(l — P)a, and a\ 
= (p + q)Xp/q — pX/q . Collisions of the particle with the 
planet occur on dashed lines in the kh plane. The dashed 


lines cannot be crossed by the kh point during the evolution 
in an MMR. 


6 STATIONARY/PERIODIC SOLUTIONS 

Stationary solutions can be found as intersections of solu¬ 


tions of the second and the third equation in Eqs. (311 at 


the universal eccentricity for one chosen semimajor axis and 
for one chosen longitude of pericenter. The semimajor axis 
will be given by a shift A from an exact resonant semima¬ 
jor axis determined by the condition a T = ap (1 — P) 1 ^ 3 
[Mi,/{Mi, + Mp)] 1//3 \p/{p + q)] 2/3 - From this definition we 
have A = a — a T . Solutions of the second equation (red 
line) and third equation (green line) in Eqs. ( |31| > are de¬ 
picted in Fig. [2] Collisions of the particle with the planet 
occur on curves depicted with dashed lines. Solved system 
comprises the Sun and the Earth as two main bodies in the 
planar CRTBP with radiation. The equations were solved 
for all dust particles with Q' pT = 1 which can be in bound 
orbits around the Sun (p < 1) and are captured in the ex¬ 
terior mean motion 6/5 orbital resonance with the Earth. 
The eccentricity of particle orbit is equal to the universal 
eccentricity and the longitude of pericenter of particle orbit 
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Figure 5. Stationary solutions found from Eqs. ( |31[ ) ■ Dependencies of the shift from the exact resonance A on 3 (left column) and 
dependencies of /3 on the resonant angular variable a (right column) are shown for 6/5, 7/6, 8/7 and 9/8 exterior MMRs in a circular 
Sun-Earth-dust problem with the PR effect and the radial solar wind. Equal color in a given row for a given resonance corresponds to 
the same sets of the stationary solutions. The semimajor axis of the solutions can be computed as a = a r + A. The stationary solutions 
in the exact resonances are marked with black circles. The stationary eccentricity is equal to the universal eccentricity. 


is equal to zero for both panels. In reality, the positions of 
the lines in Fig. [2] should not be depended on the longitude 
of pericenter (see Section [4|. Depending on used numerical 
averaging method their positions can be found to be slightly 
depended on u. The left panel has A = 0 and the right panel 
has A = 3 x 10 -6 AU. As can be seen in Fig. [ 2 ] the posi¬ 
tions of the green lines are not significantly dependent on A 
while the positions of the red lines are. Top boundaries for 
the collisions in Fig. [2] are obtained when an aphelion of the 
particle orbit is located on the planetary orbit. The value of 
ft op on the boundaries can be calculated as 


ft op = 1 - 


M* + M P 
M * 


p + q 


1 + eii, 



, (32) 


where eii m is the universal eccentricity. Solutions of the 
third equation (green lines) for large values of /3 are close 
to the collisions. At large fi strong gravitational pull from 
the planet is needed in order to compensate the influence 
caused by the PR effect and the solar wind (Eqs. 291. Very 
large values of /3 are shown in Fig.[2]only for completeness of 
the depicted solutions since a capture of the dust particle in 
an MMR at these values of /3 is complicated by strong per¬ 
turbation from the non-gravitational effects. The stationary 
solutions are located in intersections of the red lines and the 
green lines in Fig. [2] The intersections are marked with black 
circles. Two intersections exist for each panel. The top in¬ 
tersections (subscript t) are located close to collisions while 
the bottom intersections are not. 
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Figure 6. Pericenters of orbits corresponding to the stationary solutions shown in the reference frame orbiting with the Earth for dust 
particles captured in 6/5, 7/6, 8/7 and 9/8 exterior mean motion orbital resonances in a circular Sun-Earth-dust problem with radiation. 
The eccentricity of orbits is equal to the universal eccentricity. Intervals for Gj are calculated from g depicted in the right column of Fig. 
[5] (see text). /3 max is maximal /3 at which the stationary solution exists and /3 max belongs to the 6/5 resonance (/9 max — 0- 532394). 


Orbital evolution obtained from numerical solution of 
the equation of motion (Eq. |26[ ) during a stationary solu¬ 
tion is presented in Fig. [3] All evolving parameter are ob¬ 
tained from averaging over the synodic period (evolution of 
osculating orbital parameters during initial 30 years of this 
numerical solution can be seen in Fig. |3|. The stationary 
evolution is found for the stationary point Eb depicted in 
Fig. i Numerically found A of the averaged semimajor axis 
is 3.0057 x 10 -f> AU. This is in an excellent agreement with 
value 3 x 10 -6 AU predicted by presented analytical theory 
for this dust particle (see Fig. [2j. During the time of integra¬ 
tion the longitude of pericenter increased by about 4.2419°. 
At greater resolution of vertical axes in Fig. [3]it is possible to 
see oscillations in the evolutions of a, e and a. This is caused 
by finite possibilities of the numerical integration. The sta¬ 
tionary e and a determined by the numerical solution of the 
equation of motion are found to be very sensitive on a pre¬ 
cision of the numerical integration. Differences between the 


values of the stationary e and a obtained from the presented 
analytical theory and the values obtained from the numer¬ 
ical integration can be caused also by limits of averaging 
used in the determination of the analytical values. 

We searched the stationary solutions for exterior reso¬ 
nances 6/5, 7/6, 8/7 and 9/8 in a circular Sun-Earth-dust 
problem with the PR effect and the radial solar wind. We 
have chosen these exterior resonances since the orbits corre¬ 
sponding to the stationary solutions for exterior resonances 
of the first order (q = — 1) with p < 6 reach close to or be¬ 
hind the orbit of Mars for particles with small /3. We must 
note that the orbits corresponding to the stationary solu¬ 
tions for the considered exterior resonances at large /? (~ 
0.5) can reach the orbit of Venus. The stationary solutions 
are shown in Fig. [5] Dependencies of A on ^ are shown in 
the left column of the figure and dependencies of /3 on er are 
shown in the right column of the figure. Equal color on both 
sides corresponds to equal set of the solutions. We will call 
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these sets of the solutions branches. The bottom branches 
of the left-hand side plots begins at /3 — 0 and ends with the 
solution with A = 0. The stationary solutions for the exact 
resonances are marked with black circles. The top branches 
of the left-hand side plots ends with minimal /? for which 
a crossing of collision curves occur if are all stationary A 
considered. Such a crossing is located in Fig. [2] close to the 
left and the right margins of both plots. The crossing of col¬ 
lision lines corresponds to such a particle on such an orbit 
that collision occurs during one orbit before and also after 
the pericenter of the particle orbit (in Keplerian approxima¬ 
tion of the motion). 

According to the results in Section [4] the angle between 
the planet and the pericenter of the particle orbit does not 
change after the synodic period. In Fig.[6]are shown pericen- 
ters of the stationary solutions in a reference frame orbiting 
with the Earth. The number of stationary branches in this 
reference frame is multiplied by p + q. This is caused by the 
periodicity of a with the period 2n and by the identity ui = 
— q / (p + q) u holding for the particle in pericenter and the 
planet located on the x-axis (Eq. 12). In the Fig. [6] we can 
see that the pericentes create a ring inside of Earth’s orbit. 
When the particle is in pericenter, then a difference between 
angular velocities of the planet and the particle is minimal. 
The difference can be also negative for the considered ex¬ 
terior mean motion resonances. Therefore the particle can 
spend longer time in these places with respect to the planet 
( Kuchner fe Holman||2003 1. 

The stationary solutions for the considered resonances 
located on the bottom branches in the left plots in Fig.[5]can 
contribute to a cloud trailing the Earth observed by satellites 
IRAS ( |Dermott et al.|1994 i and COBE ( Reach et al.|1995} . 
These branches are located immediately below the Earth in 
Fig. The contribution from these resonances should lead 
to a trailing cloud which should occupy large interval of he¬ 
liocentric distances. This is inconsistent with observed shape 
of the cloud (see Fig. 2b in Reach et al. 1995[ ). Narrow cir¬ 
cumsolar dust ring was observed also close to the orbit of 


Venus by the STEREO mission (Jones, Bewsher & Brown 
[2013}. A special case mean motion 1/1 orbital resonance 
produces narrow circumstellar rings since the universal ec¬ 
centricity for this resonance is zero. The particles captured 
in the 1/1 resonance can stay close to the same place with 
respect to the planet and produce the trailing cloud. How¬ 
ever there is problem with capture the dust particles into the 
1/1 resonance. It is practically impossible for a dust particle 
to drift towards and get trapped in the 1/1 resonance when 
the dust particle is approaching the Sun under the PR effect 
and radial stellar wind (Liou et al. 1 995} . The capture into 
exterior resonances is possible. 


7 CONCLUSIONS 

We have investigated stationary solutions of dust orbits for 
cosmic dust particles under the action of the PR. effect and 
the radial stellar wind captured in MMRs with a planet 


in a circular orbit. The secular resonant evolution can be 
described using averaged resonant equations. The averaged 
resonant equations yield a system of equations valid for sta¬ 
tionary points. 

The stationarity of the semimajor axis and of the li- 


bration center (kh point) is most often considered dBeauge 


fe Ferraz-Mello|l993 1994 Sidlichovsky fe Nesvornj'J1994 1 . 
The validity of the conditions for the stationary point in 
astrophysical systems with rotational symmetry yields a pe¬ 
riodic motion in the reference frame orbiting with the planet 
around the star. 

We successfully solved the averaged resonant equations 
restricted by the conditions for the stationary points when 
the non-gravitational effects are the PR effect and the radial 
stellar wind. Search of the stationary points using numerical 
solution of the equation of motion was successful. Numeri¬ 
cally found shift from the exact resonance was in an excellent 
agreement with the value predicted by the analytical theory. 

In dynamical modeling of the debris disks is improb¬ 
able that the dust particle captured to the exterior MMR 
after random capture process obtains conditions leading to a 
permanent capture with the zero libration amplitude (after 
averaging over the synodic period). However, the dust par¬ 
ticles can be placed in the vicinity of such initial conditions 
(see Fig. [3}. All captures in the exterior MMRs in the planar 
CRTBP are unstable in the sense that the libration ampli¬ 


tude increases in accordance with results of Gomes (19951. 


Even small libration amplitude will be finally increased and 
the capture leads to a close encounter with the planet (in 
Fig.[l]the kh point crosses the dashed line). The periodicity 
discussed in this study belongs to exact initial conditions 
that give the zero libration amplitude. The zero libration 
amplitude does not increase ( Gomes| 19951. 

Interesting result is that the existence of periodic mo¬ 
tions in the reference frame orbiting with the planet results 
from the averaged theory. The existence of the periodical 
motions in the MMRs occurring in the planar CRTBP with 
radiation demonstrated by the existence of circular orbits for 
the 1/1 resonance ( Pastor| 2014b'| would be nicely completed 
by the existence of the periodical motions in the exterior 
resonances. The stationary solution for a given dust particle 
captured in a given exterior MMR in the planar CRTBP 
with radiation would represent a core that is a periodic mo¬ 
tion to which the eccentricity asymptotically approaches and 
around which the libration occurs. 


APPENDIX A: SECULAR TIME DERIVATIVES 
IN ROTATIONAL SYMMETRY SYSTEMS 


The Gaussian perturbation equations of celestial mechanics 
are (cf. e.g. Danby 1988 |Bate et al.||197l ) 


da 

dt 

de 

dt 


2 

ny/1 — e 2 

an 


: [oRe sin / + ax (1 + e cos /)] 
a R sin / + a T ( cos / 
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duj 

dt 


da b dn 

—- + t — 

dt dt 


+ 


e + cos / 
1 + e cos / 


\/l-i 


— ax sin / 


1 - e 
an 

— Ot 


or 


a R cos / 


2 + e cos / 
1 + e cos / 
/ cos / _ 


1 + e cos / 


sin / 2 + e cos / 
e 1 + e cos / 


(Al) 


where / is the true anomaly, a R is the radial component of 
the perturbing acceleration and ax is the transversal com¬ 
ponent of the perturbing acceleration. Let the particle be in 
some position in the Keplerian orbit which has the longitude 
of pericenter cJi. For the perturbing acceleration of a dust 
particle caused by the non-gravitational effects we have in 
Cartesian coordinates defined in the orbital plane 


api = (apn, apiy) . 


(A2) 


Orthogonal radial and transversal unit vectors of the dust 
particle in the orbit are 


e R i = [cos(/ + wi),sin(/ + wi)] , 

exi = [-sin (/+ wi), cos (/-t-wi)] . (A3) 


We can calculate a R i and axi using 
ffl R i = api • eki = apix cos (f + u>i) + a P i w sin (f + u> i) , 

axi = api • exi = — ap\ x sin (/ + wi) + api H cos (/ + un) ■ 

(A4) 

Since the problem has the rotational symmetry it is possible 
to obtain the perturbing acceleration in an orbit with dif¬ 
ferent longitude of pericenter ©2 by a rotation around the 
origin of the Cartesian coordinates by the angle cJi — C 02 - 

ap 2 x = apix cos (oji — UI 2 ) + apiy sin (wi — lu 2 ) , 

ap 2 y = — apix sin (ah — 1 I 12 ) + a P i y cos (lji — aj 2 ) ■ (A5) 


The radial component of the perturbing acceleration for the 
particle on the orbit with the longitude of pericenter ui 2 is 


OR2 = CtP2 ■ e R 2 

= [ap\x cos (wi — W 2 ) + apiy sin (uii — £ 2 )] 

X COS (/ + UO2) 

+ [—apix sin (wi — w 2 ) + apiy cos (w 1 — 0 J 2 )] 
x sin (/ + 0 J 2 ) 

= apix cos (/ + wi) + apiy sin (/ + wi) = a R i . (A6) 


Similarly for the transversal component of the perturbing 
acceleration we obtain 


ax 2 = axi ■ (A7) 

We come to the conclusion that the radial and transversal 
components of the perturbing acceleration are independent 
of the longitude of pericenter. Because the longitude of peri¬ 
center is not explicitly present on the right-hand sides of 


the Gaussian perturbation equations given by Eqs. (Al 1 , af¬ 
ter the averaging process we obtain the secular expressions 
independent of the longitude of pericenter. 
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